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Abstract. We consider a Susceptible-Infective-Recovered (SIR) model, where the mechanism for the re- 
newal of susceptibles is demographic, on a ring with next nearest neighbour interactions, and a family of 
correlated pair approximations (CPA), parametrized by a measure of the relative contributions of loops 
and open triplets of the sites involved in the infection process. We have found that the phase diagram of the 
CPA, at fixed coordination number, changes qualitatively as the relative weight of the loops increases, from 
the phase diagram of the uncorrelated pair approximation to phase diagrams typical of one-dimensional 
systems. In addition, we have performed computer simulations of the same model and shown that while 
the CPA with a constant correlation parameter cannot describe the global behaviour of the model, a rea- 
sonable description of the endemic equilibria as well as of the phase diagram may be obtained by allowing 
the parameter to depend on the demographic rate. 

PACS. 02.50.-r Probability theory, stochastic processes, and statistics - 87.23.Ge Dynamics of social 
systems - 05.70.Ln Nonequilibrium and irreversible thermodynamics 



\ 1 Introduction 

' Stochastic Susceptible-Infective- Recovered (SIR) epidemic 
[ models on lattices and networks can be mapped on to per- 
. eolation problems and are well understood laai To 
' describe disease spread and persistence in a community, 
, the model must be extended to include a mechanism for 
' renewal of susceptibles, either births or immunity waning. 

Models with immunity waning, Susceptible-Infective- 
' Recovered-Susceptible (SIRS), are based on the following 
I transitions: 

S^^I^^R^^S, (1) 

In 

meaning that any susceptible individual S can be infected 
by an infected neighbour at the infection rate /3, any 
infected individual / becomes recovered R at the recovery 
rate S, and any recovered individual R becomes susceptible 
iS* at the immunity loss rate fj.. Following customary habits, 
we shall choose time units for which 6—1. 

The SIRS model interpolates between two well known 
models, the contact process (also known as Susccptible- 
Infective-Susceptible or SIS) and the SIR model, in the 
limits /i — > oo and /i — > 0, respectively, and much is known 
about its behaviour on regular lattices, both from the 
point of view of rigorous results 0, |E Hi and of assessing 
the performance of mean field and pair approximations 
against stochastic simulations 0. 
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In particular, it is known that on hypercubic lat- 
tices of arbitrary dimension the phase diagram of ^ has 
two critical values, Pc{oo) < /3c(0), which are the critical 
rates of the two limit problems that is the contact process 
and SIR, respectively. For f3 < Pc{oo) there is disease ex- 
tinction for every y^t, while for /3c (0) < /3 there is disease 
persistence for every /i. For /3c(oo) < (3 < /3c(0) disease 
persistence occurs only for /i above a certain threshold. 
The region of disease persistence for every ^ is 'missing' 
in dimension (i=l, because in this case /3c(0) is infinite. 

In the uncorrelated pair approximation (UPA, see 
Section|2j) was applied to the SIRS model on linear and 
square lattices, and the phase diagrams computed from 
the corresponding equations of evolution were compared 
with the mean field phase diagram and with the results 
of simulations. It was shown that, by contrast with the 
mean field approximation, the UPA phase diagram agrees 
qualitatively with the simulations and the exact results 
both in d= \ and in d = 2. 

Since the UPA does not take into account the lattice 
dimensionality explicitly, it predicts identical phase dia- 
grams on lattices with the same coordination number fc, 
namely on linear ((i=l) and square (d = 2) lattices, when 
next nearest neighbours (fc = 4) are considered. However, 
in one dimension the critical infection rate /3c (0) = oo, and 
the critical line has an asymptote at = 0, while in two 
dimensions the critical line crosses the fi = Q axis at a 
finite value of (3c^ which is the result of the UPA for k = A. 

The question is then whether generalized pair approx- 
imations can account for the dependence on dimension- 
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ality and. in particular, whether they can describe phase 
diagrams with different quahtative behaviours at fixed co- 
ordination numbers. 

We have addressed this question, and more gcneraUy 
the problem of constructing suitable pair approximations 
(Section 13), m the context of a modification of model (0), 
where the mechanism of renewal of susceptibles is demog- 
raphy, rather than immunity waning. This is the natural 
scenario in the epidemiology of diseases that confer per- 
manent immunity, such as childhood infectious diseases 
1^, . For this model infection obeys the same rules as 
in immunity is permanent and all individuals, what- 
ever their state, are submitted to birth and death events 
at a rate fi. The stochastic process, which describes the 
dynamics of this system, is governed by the transitions 



(a) = 1.00 



(b) fj. = 10.0 



S^I^R, 



{S,I,R}-^S. 



(2a) 
(2b) 



In the limit /-t = 0, both models and (0) coincide with 
SIR model. In the opposite limit the dynamics of the two 
models are drastically different. While, in the limit /i = oo, 
SIRS coincides with the contact model 0, in the same 
limit the dynamics of (O is trivial: it is driven by demog- 
raphy, that keeps the entire population susceptible for any 
Pc, and thus (3c{oo) = oo. We are interested in the regime, 
where /i is smaller than the recovery rate, which is mean- 
ingful for the study of acute disease spread. Although in 
this regime the dynamics is dominated by the infection 
and recovery processes which are identical in both mod- 
els, the behaviour of 101 appears to be different, in a subtle 
way, from that of (Section [SJ. 

We have considered the demographic SIR model on 
a linear lattice with periodic boundary conditions (ring) 
and next nearest neighbour interactions, A: = 4. We con- 
structed a family of correlated pair approximations (CPA), 
parametrized by 0, a measure of the relative contributions 
of loops and open triplets of connected sites involved in 
the disease spread (Sectional). For 6 ~ Q the approxima- 
tion reduces to the standard UFA (Scction|21). The phase 
diagrams of the CFA show that as 9 increases from to 
6* (see Section|31) the CFA interpolates between the fc = 4 
UFA critical behaviour and the typical one-dimensional 
phase behaviour, with /?c(0) = oo. Finally, we have simu- 
lated the demographic SIR model JSJ on a ring, with k~A. 
The results of the simulations indicate that while the CFA 
with a constant value of 9 cannot describe the global phase 
diagram of ij^J , a reasonable description of endemic equi- 
libria as well as of the phase diagram is obtained when 
6 is allowed to depend on the demographic rate fji (Sec- 
tionOJ. This illustrates that in addition to describe the di- 
mensional crossover for lattices with coordination number 
fc = 4, the CPA can be made semi-quantitative providing 
an alternative to the stochastic simulations of individual 
based models. We conclude in Section 0| with a brief dis- 
cussion of the results. 




Fig. 1. Endemic infective probability versus infection rate at 
(a) high demographic rate, and at (b) huge demographic rate: 
the endemic infective probability is plotted from simulations 
(open circles), the MFA (long dashed lines), the UPA (dashed 
dotted lines) and the correlated model with best-fit closed form 
parameters (solid lines). The fitting procedure is based on per- 
pendicular offsets and on the assumption that the closed form 
parameter 9 depends only on the demographic rate /i. Closed 
form parameters 9 for (a) and (b) respectively: 0.50, 0.70. 

2 Mean Field and Uncorrelated Pair 
Approximations 

In this section we consider the time evolution of the de- 
mographic SIR model on regular lattices and review the 
mean-field and (standard) uncorrelated pair approxima- 
tions, setting the notation and the stage for the develop- 
ment of the more sophisticated correlated pair approxi- 
mations. 

In the demographic SIR model on networks, sites rep- 
resent individuals and bonds social links. The dynamics is 
governed by the stochastic process (jJl . Denoting by Pt {A) 
the probability for an individual to be in state A (at time 
t), Pt{AB) the probability for a lattice bond to connect 
an individual in state A to an individual in state B, the 
time evolution of the singleton probabilities Pt [A) can be 
described by the set of first order differential equations 
Hi: 



^Pt{S) 
dt 

mil) 

dt 

dPtjR) 
dt 



-/i [Pt{r) + Pt{R)]-l3Y,Pt (SI. ) , (3a) 

n 

-|3Y,Pt{SIn)-{^^ + 5)Pt{I), (3b) 

n 

-6Pt{I)- ^Pt{R), (3c) 



where the summations run over the connected neighbours. 
Clearly the set of equations l|3Jl is not closed since it in- 
volves pair probabilities without describing their time evo- 
lution. This follows from the stochastic process where 
infection H2a|) proceeds via SI contact pairs. 

As a matter of fact, the time evolution of the g-tuple 
probabilities is described by a set of first order differential 
equations expressing their time derivatives as linear com- 
binations of g-tuple and ((7-|-l)-tuple probabilities, subject 
to a normalization condition. In order to proceed, the set 
of equations must be closed, that is the (g+l)-tuple prob- 
abilities must be written in terms of q-tuple probabilities. 
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The 'art' is to use closures that capture key physical fea- 
tures of the system and are still manageable by symbolic 
or numerical-symbolic computation. The results of a par- 
ticular closure, or approximation, may then be checked 
against rigorous results and/or stochastic simulations. 

For most closures the (g + lj-tuple probabilities are 
rational functions of the g-tuple probabilities, appropri- 
ately normalized, and thus the constrained set of first 
order differential equations may be replaced by an un- 
constrained set where the time derivatives of independent 
g-tuple probabilities are expressed as rational functions 
of these g-tuple probabilities. Although the resulting sets 
of equations are easily integrable by classical numerical 
methods and admit polynomial systems as steady state 
equations, their analysis remains cumbersome even at low 
order q. 

The simplest closure is the mean field approximation 
(MFA), where the pairs (2-tuples) are assumed to be formed 
by uncorrelated singletons (1-tuples): 

Y,PtiSI,.)^kPt{S)PtiI). (4) 

n 

For the demographic SIR model the endemic equilibrium 
(steady state) is computed easily. The mean-field endemic 
infective probabilities are plotted in Figure^as a function 
of the infection rate, at two values of fi. For any value of fi, 
the MFA predicts two different steady states: at infection 
rates /3 smaller than the critical infection rate Pc there is 
disease extinction, while at infection rates /3 greater than 
the critical infection rate Pc there is disease persistence, 
i.e. infected (and recovered) individuals coexist with sus- 
ceptibles. The two regimes are separated by the mean-field 
endemic threshold that is plotted in FigureEl (dashed hne). 

We anticipate that the results of the MFA will be accu- 
rate when the demographic process of ^ dominates over 




Fig. 2. Phase diagram for the UPA: the no-coexistence phase 
and the coexistence phase are separated by the critical curve 
from simulations (open circles), the MFA (long dashed line), 
and the UPA (thick solid line). Within the coexistence phase, 
at very low demographic rates fi, the UPA predicts an oscilla- 
tory phase as shown in the inset. 



the infectious one since in this regime pairs are continually 
broken and thus the behaviour of each individual is essen- 
tially independent on that of the other ones. The infection 
process governed by Susceptible-Infective contact pairs, 
dominates in the opposite regime (^ <C 1), relevant in the 
epidemiological context. The appropriate mean field the- 
ory is then the uncorrelated pair approximation (UPA). 

The UPA is for pairs what the MFA is for singletons. 
In the UPA triplets (3-tuples) are assumed to be formed 
by uncorrelated pairs: 

The UPA is expected to outperform the MFA but, in 
general, its solution is not known in closed form. For the 
demographic SIR model the calculation of the phase dia- 
gram and the stability analysis is still tractable by sym- 
bolic computation. For lattices with coordination number 
A; = 4, the phase diagram is plotted in Figure 13 It is clear 
that the UPA is quantitatively superior to the MFA when 
compared with the results of simulations (open circles). 
Both the MF and the UP approximations of the fc = 4 
demographic SIR model predict a finite critical infection 
rate at /i = 0, while the simulations indicate that f3c will 
diverge as fi tends to 0. 

However, the SIRS and the demographic SIR models 
are different at low (but finite) demographic rates fi. In 
the demographic SIR model the mechanism for the re- 
newal of susceptibles is totally random by contrast to the 
mechanism of the SIRS model. In our model susceptibles 
are born anywhere on the lattice while in the SIRS model 
only previously infected sites loose immunity. 

We note that the randomizing effect of the demographic 
SIR mechanism for the renewal of susceptibles is reminis- 
cent of the randomizing effect of shortcuts in small-world 
networks of the Watts and Strogatz type 0, ^| where 
correlations are destroyed and an effective mixing of the 
population is achieved, with (drastic) consequences on the 
phase diagram. 

Finally, it is worth noticing that the UPA also predicts 
the existence of an oscillatory phase within the survival 
or coexistence phase {i.e., to the right of /3c (0)), for small 
values of /i (Figure EJ. The same is true for the UPA of 
process Q on the square lattice. This behaviour will be 
difficult to identify in stochastic simulations, since it may 
be blurred by large fluctuations and stochastic extinctions. 

3 Correlated pair approximations 

In order to construct more realistic pair approximations, 
we have investigated closure procedures inspired by the 
geometrical structure of the lattice. 

Within this perspective and as far as social triplets 
are concerned, the ring of degree fc = 4 and the triangular 
lattice (A: = 6) are propitious networks since their nearest- 
neighbour triplets split into two distinct classes: 'chain- 
like' (open) and 'loop-like' (closed) triplets. A very naive 
idea is to take into account the two classes of triplets and 
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Fig. 3. Isoparametric phase diagrams for the correlated model; 
the no-coexistence and coexistence phases are separated by the 
critical curve from simulations (open circles), the MFA (long 
dashed line), the UPA (dashed dotted line) and the correlated 
model for different 6* (solid lines). For 6** « 0.3807 (bold solid 
line) the critical infection rate jSc tends asymptotically to in- 
finity when the demographic rate vanishes. Closed form pa- 
rameters 9 from left to right: \, 9* , 
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to use the probability 6 and 1 — of finding respectively a 
'loop-like' triplet and a 'chain-like' triplet as a parameter 
to be fitted to simulation results. Thus triplets are as- 
sumed to be formed either of uncorrelated (chained) pairs 
or of correlated (looped) pairs : 



' (h ^^\\(^ - ft\ MSA)Pt(SI) r, Pt(AI)Pt(SA)Pt 
^) [l^ p^j^j + U p^(A)Pt(S)Pt{ 

if A e {S, R}, 
(k - 1) PtiSI) - En [Pt{SSI„)+PtiRSQ] 
if A = /. 



(SI) 

(I) 



(6) 



The demographic SIR version of the CPA ((HJ is ame- 
nable by cumbersome numerical-symbolic computation al- 
though some interesting results may be obtained by sym- 
bolic computation. The phase diagrams are shown in Fig- 
ureH We find that, as 9 increases from to 6'*w0.3807\ 
keeping k = A fixed, the CPA phase diagrams interpolate 
between the UPA behaviour and typical one-dimensional 
phase diagrams with /3c{Q) = oo- At 6* the critical in- 
fection rate Pc tends asymptotically to infinity as the de- 
mographic rate ^ vanishes. Inspection of Figure 13 also 
shows that the closed form parameter 9 cannot be con- 
stant if a quantitative description of the global phase dia- 
gram is required. If we allow 6 to depend on ^, reasonable 
descriptions of the endemic equilibria (Figure ^ and of 
the global phase diagram (Figure^ are obtained. For the 



^ the real solution of the cubic equation 276'^ - 189"^ + 879 
32 = 0. 



SIRS model on the square lattice, a CPA obtained 
by fitting 6* to /?c(0) will improve the results of the UPA 
used in Q to describe the behaviour of the system at low 
values of ^. 



4 Discussion 

We have proposed a simple CPA that was shown to pro- 
vide a reasonable approximation to the behaviour of stochas- 
tic models that are relevant in epidemiology — the agree- 
ment against simulation data being far better than MFA 
and UPA with a suitable choice of the parameters. The 
resulting equations of evolution may be used to approxi- 
mate phase diagrams, as well as steady state and dynam- 
ical behaviours of the associated stochastic models. The 
CPA takes into account some of the effects of the local 
lattice structure and yields a clear alternative to heavy 
stochastic simulations. 

One of the directions of future work includes the de- 
velopment of CPAs, along the lines of the present work, 
to account for the local (lattice like) structure of a class of 
complex networks, such as the Watts and Strogatz small- 
world networks, that have been shown to be relevant in 
epidemiological contexts • 
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